rm(list=ls())
setwd('~/Dropbox/Farm Subsidies and Committees/analyses/local_reg/')

library(ggplot2)
library(dplyr)
library(magrittr)
library(stringr)
library(haven)
library(statar)
library(gridExtra)


for (switch in c('all', 'lfp', 'dcp')) {
    df = read_dta(paste0('local_reg_', switch, '.dta'))
    base_dvs = c('voted', 'ran', 'won', 'woncondran')
    dvs = c(paste(base_dvs, 'avg', sep = '_'), paste(base_dvs, 'max', sep = '_'))
    
    # get bin data
    for (dv in dvs) {
        tmp_gg = ggplot(data = df, mapping = aes_string(x = 'amt_avg_wsz/1000', y = dv)) + stat_binmean(n = 100)
    
        tmp_bins = layer_data(tmp_gg) %>%
            group_by(x) %>%
            summarise(y = mean(y))
        
        assign(paste('bins', dv, sep = '_'), tmp_bins)
    }
    
    
    # plots
    for (dv in dvs) {
        ## title
        tmp_title = ''
        if (str_detect(dv, '^voted')) {
            tmp_ylab = 'Voted'
        } else if (str_detect(dv, '^ran')) {
            tmp_ylab = 'Ran'
        } else if (str_detect(dv, '^won_')) {
            tmp_ylab = 'Won'
        }
            
        if (str_detect(dv, '_avg$')) {
            tmp_title = paste(tmp_ylab, '(mean across years within ID)')
        } else if (str_detect(dv, '_max$')) {
            tmp_title = paste(tmp_ylab, '(max across years within ID)')
        }
        
        ## plot (thousands)
        tmp_bins = eval(as.name(paste('bins', dv, sep = '_')))
        
        tmp_gg = ggplot(data = df, 
                        mapping = aes_string(x = 'amt_avg_wsz/1000', y = dv)) + 
            stat_smooth(data = df, 
                        mapping = aes_string(x = 'amt_avg_wsz/1000', y = dv), 
                        method = 'gam', 
                        formula = y ~ s(x, bs = 'cs'), 
                        color = 'red', 
                        linetype = 'dashed', 
                        size = 0.3) +
            geom_point(data = tmp_bins, 
                       mapping = aes(x = x, y = y)) +
            labs(x = 'Mean amount disbursed (thousands)', 
                 y = paste(tmp_ylab, '(percent)')) + 
            theme_minimal() + 
            theme(plot.title = element_text(hjust = 0.5),
                  plot.caption = element_text(hjust = 0))
        
        assign(paste('plot', dv, sep = '_'), tmp_gg)
    }
    
    
    note = '\nNote: Points are means of DV plotted against means of mean payment amount within 50 bins. Smoothed over entire data using a generalized additive model. Average payment amount is winsorized at the 1st and 99th percentiles.'
    note %<>% strwrap(160) %>% paste0(sep = '', collapse = '\n')
    
    ggsave(paste0('grid_avg_', switch, '.pdf'), 
           arrangeGrob(plot_voted_avg, plot_ran_avg, plot_won_avg, nrow = 1),
           height = 3, width = 9)
    
    ggsave(paste0('grid_max_', switch, '.pdf'), 
           arrangeGrob(plot_voted_max, plot_ran_max, plot_won_max, nrow = 1),
           height = 3, width = 9)
}


